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Case  History 

Imaging  a  hydrate-related  cold  vent  offshore  Vancouver  Island  from 
deep-towed  multichannel  seismic  data 

Tao  George  D.  Spence^,  Warren  T.  Wood^,  Michael  Riedel'^,  and  Roy  D.  Hyndman®'^ 


ABSTRACT 

The  Bullseye  vent,  an  approximately  500-m-diameter  deep- 
sea,  hydrate-related  cold  vent  on  the  midslope  offshore  Vancou¬ 
ver  Island,  was  imaged  in  a  high-resolution  multichannel  survey 
by  the  Deep-towed  Acoustics  and  Geophysics  System  (DTAGS) 
The  structure  was  drilled  by  the  Integrated  Ocean  Drilling  Pro¬ 
gram  at  site  U1 328.  Towed  about  300  m  above  the  seafloor,  the 
high-frequency  (220-820  Hz)  DTAGS  system  provides  a  high 
vertical  and  lateral  resolution  image.  The  major  problems  in  im¬ 
aging  with  DTAGS  data  are  nonlinear  variations  of  the  source 
depths  and  receiver  locations.  The  high-frequency,  short- wave¬ 
length  data  require  very  accurate  positioning  of  source  and  re¬ 
ceivers  for  stacking  and  velocity  analyses.  New  routines  were  de¬ 
veloped  for  optimal  processing,  including  receiver  cable  geome¬ 
try  estimation  from  node  depths,  direct  arrivals  and  sea-surface 
reflections  using  a  genetic  algorithm  inversion  method,  and 
acoustic  image  stitching  based  on  relative  source  positioning  by 


crosscorrelating  redundant  data  between  two  adjacent  shots. 
Semblance  seismic  velocity  analysis  was  applied  to  common-re- 
flection-point  bins  of  the  corrected  data.  The  processed  images 
resolve  many  subvertical  zones  of  low  seismic  reflectivity  and 
fine  details  of  subseafloor  sediment  structure.  At  the  Bullseye 
vent,  where  a  40- m- thick  near-surface  massive  hydrate  layer  was 
drilled  at  U1 328,  the  images  resolve  the  upper  part  of  the  layer  as 
a  dipping  high-reflectivity  zone,  likely  corresponding  to  a  frac¬ 
ture  zone.  Velocity  analyses  were  not  possible  in  the  vent  struc¬ 
ture  but  were  obtained  1 80-270  m  to  either  side.  Normal  veloci¬ 
ties  are  in  the  upper  50  m,  but  over  the  interval  from  50  to  1 00  m 
below  the  seafloor  at  the  northeast  side,  the  velocities  are  higher 
than  the  average  normal  slope  sediment  velocity  of  approximate¬ 
ly  1590  m/s.  These  high  velocities  are  probably  related  to  the 
high  reflectivity  zone  and  to  the  bottom  portion  of  the  massive 
hydrate  detected  by  resistivity  measurements  in  the  upper  40  m 
atU1328. 


INTRODUCTION 

At  the  northern  Cascadia  margin,  the  incoming  subducting  sedi¬ 
mentary  cover  is  almost  completely  scraped  off  the  oceanic  plate  at 
the  deformation  front  along  thrusts  that  extend  down  near  the  top  of 
the  ocean  crust,  forming  a  large  accretionary  prism  (Davis  and 
Hyndman,  1989).  Within  the  accretionary  prism  off  Vancouver 
Island,  a  cold  vent  field  in  the  midslope  region  has  been  studied  in¬ 
tensively  since  1999,  using  many  geophysical,  geological,  and 
geochemical  techniques.  The  most  prominent  vent,  the  Bullseye 
structure,  has  a  diameter  of  approximately  500  m.  It  was  drilled 
by  the  Integrated  Ocean  Drilling  Program  (lODP)  at  site  U1328 


(Riedel  et  al.,  2006a),  as  shown  in  Figures  1  and  2,  where  high  con¬ 
centrations  of  gas  hydrate  were  detected  on  downhole  logs  and  re¬ 
covered  at  depths  of  about  5-40  m  below  the  seafloor  (mbsf). 

The  cold  vent  field  is  located  on  a  bathymetric  bench  at  1 200- 
1 400  m  depth.  lODP  site  U 1 328  and  nearby  Ocean  Drilling  Program 
(ODP)/IODP  sites  889/Ul  327  were  located  near  a  topographic  high 
on  this  bench  (Figure  1).  All  sediments  consist  of  silty  clays  and 
clayey  silts  interbedded  with  fine  sand  turbidities  (Westbrook  ct  al., 
1994).  Deformed,  compacted,  and  partially  cemented  Cascadia  Ba¬ 
sin  sediments  are  found  below  128  mbsf  at  site  U1327.  Little-de¬ 
formed  slope  sediments  are  found  above  128  mbsf  at  U1327  and 
above  300  m  at  U1 328  (Riedel  et  al. ,  2006a). 
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The  Cascadia  accretionary  margin  contains  one  of  the  most  com¬ 
prehensively  studied  occurrences  of  marine  gas  hydrate.  Offshore 
Vancouver  Island,  a  hydrate  bottom-stimulating  reflection  (BSR) 
occurs  over  much  of  a  30-km-wide  band  parallel  to  the  coast  beneath 
the  midcontinental  slope  (Figure  la).  As  one  end  member  for  hy¬ 
drate  formation,  methane  is  removed  from  regional  diffusive  up¬ 
ward  fluid  flow  driven  by  tectonic  compression  (Hyndman  and 
Davis,  1992;  Haacke  et  al.,  2008b).  In  these  models,  the  regional 
concentration  of  gas  hydrate  gradually  decreases  upward  from  a 
maximum  near  the  base  of  the  gas  hydrate  stability  zone  (GHSZ). 
However,  recent  results  from  expedition  311  show  that  the  hydrate 
distribution  is  very  irregular,  and  at  two  sites  the  largest  concentra¬ 
tions  of  gas  hydrate  (>80%  of  the  pore  volume  in  sections  several 
tens  of  meters  thick)  occur  at  shallow  depths  of  less  than  100  mbsf 
(sites  U1 326  and  U1327)  (Riedel  et  al.,  2006a). 

Focused  flow  within  cold  vents  represents  the  other  end-member 
mechanism  for  hydrate  formation  where  the  local  hydrate  concen¬ 
trations  may  be  very  high.  It  is  unknown  how  imp>ortant  these  cold 
vents  are  in  the  total  fluid-carbon-gas  budget  in  an  accretionary 
prism  (Riedel  et  al.,  2006a).  There  are  many  indicators  of  methane 
venting  at  the  seafloor:  gas-bubble  plumes  in  the  water  column 
(Suess  et  al.,  2001 ;  Vidalie,  2007;  E.  Willoughby,  personal  commu¬ 
nication,  2007),  chemosynthetic  communities  (Lewis  and  Co¬ 
chrane,  1 990;  Riedel  et  al.,  2006b),  massive  carbonate  crusts  (Ritger 
et  al.,  1987),  near-seafloor  gas  hydrate  (MacDonald  et  al.,  1994; 
Suess  et  al.,  2001;  Chapman  et  al.,  2004),  and  seabed  pockmarks 
(Hovland  et  al.,  1993;  Vogt  et  al.,  1999).  Seismically,  other  indica¬ 
tors  are  seismic  wipeout  zones,  or  subvertical  zones  of  reduced  re¬ 


flectivity  (Riedel  et  al.,  2002),  shallow/elevated  BSRs  (Wood  et  al., 
2004;  Trehu  et  al.,  2003),  and  local  velocity  pull-up/pull-down  fea¬ 
tures  (Scholl  et  al.,  1 978;  Lee  et  al.,  2005). 

Different  mechanisms  have  been  proposed  for  the  origin  of  the 
seismic  wipeout  and  for  the  nature  of  the  fluid  venting.  The  mecha¬ 
nism  must  explain  how  gas  or  gas-saturated  fluid  reaches  the  seaf¬ 
loor  without  forming  hydrate  as  it  enters  the  stability  field  250  m  be¬ 
low  the  seafloor.  Wood  et  al.  (2002)  propose  that  warm  rising  fluids 
significantly  perturb  the  base  of  gas  hydrate  stability  upward  to  near 
the  seafloor,  and  the  entire  blank  zone  is  from  free  gas  beneath  the  el¬ 
evated  BSR.  However,  heat-flow  measurements  have  not  shown  a 
current  thermal  anomaly  over  the  Bullseye  vent  structure,  and 
Riedel  et  al.  (2002)  suggest  that  much  of  the  wipeout  is  caused  by 
near- surface  carbonate  or  by  the  massive  hydrate  layer  at  the  sur¬ 
face;  additional  blanking  may  result  from  scattered  accumulations 
of  hydrate  in  the  subsurface.  Zuhlsdorff  and  Spiess  (2004)  argue  that 
the  blanking  zones  below  vents  are  produced  by  free  gas  bubbles 
within  natural  hydraulic  fractures.  Liu  and  Flemings  (2007)  postu¬ 
late  that  the  base  of  hydrate  stability  moves  upward  in  the  vent  struc¬ 
tures  because  of  the  progressively  higher  salinities  generated  by  the 
salt  that  moves  into  the  residual  pore  water  as  hydrate  is  formed. 
Haacke  et  al.  (2008a)  show  how  this  model  can  be  applied  to  struc¬ 
tures  off  Korea. 

We  report  on  a  deep-towed  high-resolution  multichannel  seismic 
survey  (DTAGS2),  carried  out  in  2002  over  this  cold  vent  field  near 
site  U1328  (Figure  1).  The  project  was  a  collaboration  among  the 
University  of  Victoria,  the  Geological  Survey  of  Canada,  and  the 
U.S.  Naval  Research  Laboratory  (NRL,  Wash¬ 
ington,  D.C.,  and  Stennis  Space  Center,  Missis¬ 
sippi,  U.S. A.).  The  data  were  obtained  using  the 
high-resolution  NRL  deep-towed  acoustic/geo¬ 
physical  system  (DTAGS)  on  the  CCGS  John  R 
Tally.  The  primary  goals  were  to  obtain  high-res- 
olution  structural  details  at  a  cold  vent  field,  in¬ 
cluding  the  extent  of  associated  methane  gas  hy¬ 
drate  and  underlying  free  gas,  and  to  determine 
the  associated  seismic  velocities  that  could  be 
used  to  estimate  gas  hydrate  concentrations. 

DTAGS  SYSTEM  AND  DATA 
ACQUISITION 

The  DTAGS  system  consists  of  a  Helmholtz 
resonator  source  and  a  solid-core  hydrophone 
streamer  with  two  subarrays  containing  24  omni¬ 
directional  hydrophones  (Table  1)  (Gettrustetal., 
2004).  The  advantages  of  DTAGS  over  sea-sur¬ 
face  systems  are  that  it  provides  a  source  that  is 
highly  reproducible  because  it  is  solid  state  and 
there  is  no  cavitation  and  high  vertical  resolu¬ 
tion  resulting  from  the  high  frequencies  (220- 
820  Hz).  The  wavelengths  are  less  than  6  m, 
much  shorter  than  the  25-50  m  associated  with 
typical  conventional  multichannel  seismic  data 
with  frequencies  of  30-60  Hz.  Another  advan¬ 
tage  is  greater  sampling  in  wavenumber  space  us¬ 
ing  shorter  hydrophone  arrays  with  small  receiver 
spacing.  In  addition,  a  smaller  Fresnel  zone  re¬ 
sulting  from  the  deep-tow  geometry  provides  a 
much  shorter  distance  from  the  source  to  the  re- 


Figure  1 .  (a)  General  tectonic  setting  of  northern  Cascadia  margin.  The  dashed  box  indi¬ 
cates  the  study  area  of  the  detailed  survey  location  in  the  region  of  the  Bullseye  vent  field 
on  the  midslope  region  offshore  Vancouver  Island,  (b)  Location  of  2002  DTAGS 2  lines. 
Teledyne  2000  single-channel  3D  grid  (blue  box)  and  COAMS  1999  multichannel  grid 
(green  and  gray  grid)  around  the  Bullseye  vent  (yellow  area  with  black  outline).  The  four 
seismic  blanking  zones  are  also  labeled  with  BZl  -4  as  well  as  on  Figure  14.  The  dashed 
cyan  box  outlines  the  area  showing  details  in  Figure  2.  The  bathymetry  contour  interval  is 
10  m,  using  the  data  from  a  Sonne  cruise  SO  149  in  1999  (Zuhlsdorff  etal.,  2002). 
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flectors.  Furthermore,  the  deep- towed  system  records  wider-angle 
reflection  data  and  may  provide  useful  information  for  moveout  ve¬ 
locity  and  amplitude  variation  with  offset  (AVO)  analyses  (Gettrust 
et  al.,  1 999, 2004).  The  major  disadvantage  is  that  the  high-frequen¬ 
cy,  short-wavelength  data  require  very  accurate  positioning  of 
source  and  receivers  for  stacking  and  velocity  analyses.  The  special 
processing  required  to  establish  the  position  of  the  source  and  each 
receiver  sensor  is  a  primary  focus  of  this  study. 

During  the  2002  DTAGS2  experiment,  1 1  lines  crossed  the  Bulls- 
eye  vent  field  (Figure  lb):  seven  southwest-northeast  lines  plus  four 
southeast-northwest  lines  in  the  direction.  The  source  was  fired  ev¬ 
ery  30  s,  or  ~30  m  at  a  nominal  ship  speed  of  2  knots.  Using  a 
24-bit  Sigma- Delta  converter  at  each  hydrophone,  data  were  digi¬ 
tized  at  a  sampling  interval  of  0.5  ms  with  a  recording  length  of 
1 .9  s,  including  a  250-ms  delay  (Spence  et  al.,  2002). 

Four  engineering  nodes  were  on  the  streamer:  node  1  was  on  the 
source  at  0  m,  and  the  other  three  were  on  hydrophone  channel  28  at 
138  m,  channel  38  at  288  m,  and  channel  48  (the  last  channel)  at 
438  m.  These  engineering  nodes  recorded  pressure  (depth)  and  tem¬ 
perature,  as  well  as  the  seismic  signal.  The  pressure  gauges  at  the 
four  nodes  yielded  several  pressure  readings  during  the  30  s  be¬ 
tween  shots.  These  pressures  were  averaged  over  the  shot  interval, 
converted  to  depth  (including  a  local  density  correction  based  on  a 
temperature  measurement  in  the  node)  and  placed  in  the  SEG  Y  trace 
header.  Normally  these  depths  are  more  accurate  than  the  values  cal¬ 
culated  from  traveltime  and  estimated  water  velocity.  Figure  3 
shows  the  depth  variations  of  the  nodes  during  survey  line  DT09. 

Although  the  navigation  of  the  ship  (differential  GPS,  or  DGPS) 
provided  positions  with  antenna  position  uncertainties  close  to  1  m. 


Figure  2.  Map  view  of  the  Bullseye  region,  associated  with  seismic 
blanking  zone  BZl  (light  gray);  massive  hydrate  was  recovered  in 
shallow  piston  cores  within  the  dark  gray  region  (pockmark).  The 
2002  DTAGS2  survey  lines  are  indicated  by  solid  black  lines.  The 
thick  box  outlines  the  area  of  the  2000  Teledyne  single-channel 
pseudo-3D  grid.  The  thick  gray  lines  indicate  the  seismic  sections 
showing  hydrate  cap  details  in  Figures  9-12.  The  locations  of  two 
velocity  analyses,  S3 1 7  and  S377,  are  labeled  with  black  stars.  Tran¬ 
sponder  deployment  (T1T4)  during  this  experiment  and  lODP  site 
U 1 328  are  also  labeled. 


the  DTAGS  was  towed  behind  the  ship  at  horizontal  distances  of  up 
to  1  km,  and  its  horizontal  distance  could  vary  hundreds  of  meters 
from  the  track  of  the  ship  within  a  survey.  Therefore,  navigation  was 
also  determined  acoustically  with  a  transponder  net  deployed  in 
200 1  in  preparation  for  the  DTAGS2  experiment  in  2002. 

A  four-element  long-baseline  (LBL)  acoustic  seafloor  array  was 
arranged  on  four  corners  of  a  square  with  sides  about  1  km  long 
(Figure  2).  A  fifth  relay  transponder,  or  rover,  was  attached  to  the 
DTAGS  source.  The  DGPS  from  the  ship  and  the  known  locations  of 
transponders  on  the  seafloor  were  used  via  acoustic  traveltimes  to 
determine  the  position  of  the  relay  transponder.  During  operation, 
navigation  using  all  four  transponders  was  possible  within  the  LBL 
array,  usually  out  to  100-200  m  beyond  the  array  boundaries  but 
sometimes  up  to  1  km  outside  the  array.  WINFROG  commercial 
navigation  software  calculated  real-time  positions  of  the  ship  and 
rover.  This  allowed  a  real-time  visual  assessment  of  the  positions  of 
the  ship  and  rover  (Spence  et  al.,  2002). 

Unfortunately,  many  irregular  jumps  in  the  rover  navigation, 
mainly  from  a  large  number  of  repeated  readings  and  acoustic  inter¬ 
ference  from  strong  noise  sources,  made  it  very  difficult  to  position 
the  DTAGS  accurately  and  limited  the  value  of  using  the  acoustic 


Table  1.  DTAGS  configuration  during  the  2002  DTAGS2 
cruise. 


Characteristic 

Details 

Acoustic  array 

24  channels,  group  spacing  of  3  m 

Minimum  offset:  9  m 

Maximum  offset:  78  m 

First  two  channels  dead 

Geophysical  array 

24  channels,  group  spacing  of  15  m 
Minimum  offset:  93  m 

Maximum  offset:  438  m 

Source 

Helmholtz  resonator  generating  250-ms-long 
sweep  (220-820  Hz) 

Source  strength 

—200  dB/1  Pascal  @  1  m,  sufficient  to 
penetrate  upper  400  m  of  sediments 

Shot  number 

Figure  3.  The  depth  variations  of  the  four  nodes  on  DTAGS  streamer 
along  survey  line  DT09  (Figures  lb  and  2),  showing  an  example  of 
nonlinear  variation  of  array  geometry. 
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transponder  net.  An  effort  was  made  to  correct  the  rover  navigation 
manually  in  postcruise  processing  by  eliminating  duplicate  readings 
and  acoustic  positions  that  clearly  were  in  error  and  interpolating  be¬ 
tween  reliable  positions.  The  corrected  rover  navigation  still  was  not 
accurate  enough  for  precisely  stitching  DTAGS  shot  images  as  origi¬ 
nally  designed,  but  it  provided  more  reliable  absolute  rover  3D  posi¬ 
tioning  than  navigating  the  ship  alone. 

DATA  PROCESSING 

Source  crosscorrelation  and  spectral  balancing 

The  vibroseis-like  source  of  DTAGS  is  a  250-ms-long  sweep  fre¬ 
quency-modulated  sinusoid  that  is  tapered  on  both  ends.  The  initial 
step  in  processing  the  DTAGS  data  was  to  crosscorrelate  the  field  re¬ 
cordings  with  the  source  signal.  Although  deconvolution  (frequen¬ 
cy-domain  division)  gives  abetter  spectrum  than  correlation  in  theo¬ 
ry,  correlation  (frequency  domain  multiplication)  is  in  practice  more 
stable  because  division  by  values  near  zero  will  produce  unreason¬ 
ably  large  amplitudes  at  the  corresponding  frequencies.  An  approxi¬ 
mate  pilot  sweep  was  obtained  from  the  direct  wave  recorded  at  the 
near-offset  hydrophone,  about  1 5  m  from  the  source;  previous  stud¬ 
ies  by  Chapman  et  al.  (2002)  show  that  this  approach  gives  the  clear¬ 
est  seafloor  and  deeper  reflectors. 

For  the  broadband  signals  in  the  DTAGS  data,  spectral  balancing 
achieved  better  image  quality  than  the  popular  Butterwonh  filter. 
Spectral  balancing  is  a  zero-phase  deconvolution,  effected  by  spec¬ 
tral  equalization  in  the  selected  frequency  domain;  we  assume  no 
frequency  dispersion  or  associated  phase  changes.  In  practice,  the 
2002  DTAGS  data  show  improved  image  quality  after  applying 
spectral  balancing  with  selected  band-pass  frequency  limits  of 
220-250  to  510-545  Hz:  some  seismic  events  were  enhanced  by 
compressing  the  wavelet  slightly,  and  a  relatively  clean  image  was 
achieved  by  rejecting  some  of  the  higher  frequencies  in  the  spec¬ 
trum. 


array  geometry  by  using  SSR  times  at  six  receiver  channels  and  in¬ 
terpolate  times  using  a  cubic  spline  for  the  remainder  of  the  chan¬ 
nels.  They  then  correct  traces  to  a  single  constant  datum  for  the 
whole  survey  line;  this  can  produce  large  moveout  velocity  errors  if 
the  depth  corrections  are  several  tens  of  meters. 

We  extend  the  work  of  Walia  and  Hannay  ( 1 999)  to  achieve  opti¬ 
mal  cable  geometry  and  more  accurate  time  corrections.  Our  method 
uses  a  genetic  algorithm  inversion  to  find  positions  of  the  sources 
and  receivers,  using  known  node  depths  and  manually  picked  times 
for  the  direct  arrivals  and  the  reflections  from  the  sea  surface.  Our 
method  also  uses  a  floating  datum  system  to  reduce  errors  during 
time  correction. 

We  determined  the  cable  geometry  in  the  2D  vertical  plane  be¬ 
cause  no  information  on  cable  feather  angles  was  available.  As  the 
first  step,  the  SSR  and  direct-arrival  times  were  picked  on  the  trace 
envelope  instead  of  the  raw  waveform  because  the  crosscorrelated 
wavelet  was  inconsistent  from  trace  to  trace  for  the  DTAGS  high- 
frequency  data,  caused  by  the  water-wave  roughness  of  the  sea  sur¬ 
face.  We  picked  arrival  times  at  the  first  zero  crossing  rather  than  the 
envelope  maximum,  which  shifted  with  offset  and  shot  depth  be¬ 
cause  of  phase  changes.  The  picked  times  were  then  shifted  about 
7  ms  downward  to  the  approximate  locations  of  envelope  peaks. 
Note  that  the  direct  arrival  wavelet  is  present  only  in  the  last  24  re¬ 
ceivers  at  the  farthest  offset  as  a  result  of  the  250  ms  delay  in  record¬ 
ing. 

The  second  processing  step  was  to  calculate  the  largest  horizontal 
offset  and  interpolation  velocity  for  each  shot.  Because  one  of  the 
engineering  nodes  was  attached  to  the  last  receiver,  it  was  a  unique 
point  on  the  cable  at  which  the  depth  reading,  SSR  time,  and  direct- 
arrival  time  could  all  be  obtained.  Thus,  the  horizontal  far  offset  of 
the  cable  (jc')  can  be  determined,  as  well  as  the  nominal  water  veloci¬ 
ty  (Vh,  assumed  to  be  constant  in  the  water  mass  above  the  array)  to 
be  used  for  interpolating  SSR  times  to  all  receivers  for  each  shot 
gather,  using  the  following  equations: 


Array  geometry  estimation 

Although  the  near-bottom  configuration  and  wide  high-frequency 
bandwidth  of  DTAGS  allow  much  better  resolution  of  the  subsea¬ 
floor  structures  in  the  upper  400  m,  the  system  has  three  special 
problems  for  data  imaging  and  velocity  analysis:  (1)  shot -to-shot 
variations  in  source  position  are  poorly  constrained;  (2)  difficulty 
maintaining  a  horizontal  array  shape  during  towing,  so  receiver  po¬ 
sitions  for  each  shot  are  poorly  constrained;  (3)  short  26-m  wave¬ 
lengths  that  require  very  precise  positioning  and  time  alignment  to 
determine  moveout  velocity.  Significant  nonlinear  variations  in  ar¬ 
ray  geometry,  particularly  in  the  longer  geophysical  array  (Table  1), 
occurred  throughout  most  of  the  DTAGS  experiment  (Figure  3). 
Compared  to  an  optimal  horizontal  system,  these  variations  caused 
considerable  reflection-time  errors  and  required  accurate  correction 
before  conventional  seismic  processing  could  be  applied. 

In  earlier  DTAGS  experiments  (Rowe  and  Gettrust,  1993),  depth 
readings  from  the  few  sensors  located  along  the  streamer  were  used 
to  interpolate  the  approximate  array  geometry  linearly.  However, 
these  interpolated  values  were  not  accurate  enough  to  calculate  satis¬ 
factory  time  corrections  for  the  nonlinear  cable  geometry.  The 
DTAGS  seismic  signal  normally  included  strong  sea-surface  reflec¬ 
tions  (SSR),  which  provide  more  accurate  two-way  time  and  conse¬ 
quently  depth  for  each  hydrophone  group,  assuming  constant  water 
velocity  during  a  shot  interval.  Walia  and  Hannay  (1999)  determine 


x'  =  -  idr  -  df,  (1) 


V„  = 


4dA 


^  SSR  ^  direcl 


(2) 


where  Tssr  is  the  two-way  time  of  the  sea-surface  reflection,  is 
the  one-way  time  of  direct  arrival,  and  and  dr  are  the  depth  of  the 
source  and  the  receiver,  respectively,  calculated  from  temperature- 
compensated  pressure  measured  at  the  far-offset  node. 

Assuming  that  the  depth  readings  are  accurate  and  the  cable  ge¬ 
ometry  is  smooth  to  first  order  at  locations  of  the  engineering  nodes, 
the  positions  of  traces  other  than  the  last  one  can  be  predicted  by  us¬ 
ing  their  SSR  times  to  invert  the  whole  cable  geometry  to  a  polyno¬ 
mial  spline  shape.  For  any  shot  gather,  once  the  interpolation  veloci¬ 
ty  and  farthest  horizontal  offset  JC48  at  channel  48  of  node  4  are  de¬ 
termined  by  equations  1  and  2,  a  genetic  algorithm  (GA)  method 
provided  by  MATLAB  (The  Math  Works,  Inc.)  is  used  to  determine 
the  horizontal  offsets  for  the  other  two  middle  nodes  (jcjs  of  channel 
28  and  jC3g  of  channel  38)  and  the  vertical  depths  at  three  positions: 
(diAo  dii,  and  J43  at  channels  24, 33,  and  43).  The  optimization  func¬ 
tion  of  the  GA  method  is  to  minimize  a  fitness  value  defined  as  the 
sum  Tf  of  residual  SSR  times  between  the  calculated  Tssr  and  the 
picked  Tp^k  times  for  all  channels  (channel  1 48)  within  a  shot  gather: 
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where  is  the  horizontal  offset  of  channel  /  from  the  source,  is  the 
vertical  depth  of  channel  /  from  sea  surface,  and  is  the  vertical 
depth  of  the  source  from  sea  surface. 

For  each  generation,  the  GA  method  produces  a  population  of  50 
individuals  from  the  previous  generation.  Each  individual  is  charac¬ 
terized  by  the  five  variable  values;  the  horizontal  offsets  of  the  re¬ 
maining  receivers  are  linearly  interpolated  between  the  four  nodes, 
and  then  the  depths  are  calculated  using  a  polynomial  fitting  curve. 
Order  three  and  order  five  polynomials  are  tested.  A  new  fitness  val¬ 
ue  7)  (Equation  3)  is  calculated  from  the  new  SSR  times  for  each  in¬ 
dividual,  and  the  overall  population  is  assessed  to  make  the  next  gen¬ 
eration. 

The  GA  used  for  DTAGS  cable  inversion  continues  to  produce 
generations  until  there  is  no  further  improvement  on  the  fitness 


Figure  4.  Comparison  of  the  GA  inversion  result  and  those  from  a 
linear  cable  assumption  for  shot  2230  of  survey  line  DTIO.  (a)  The 
G  A  inverted  cable  (solid  line)  has  a  more  naturally  curving  geometry 
than  the  linear  algorithm  (dashed  line  with  dots),  which  has  straight 
lines  between  cable  nodes,  (b)  Recalculated  SSR  times  from  invert¬ 
ed  cable  geometry  (solid  line)  have  much  smaller  misfit  than  the  sim¬ 
ple  linear  cable  assumption  (dashed  line)  when  compared  to  the  ac¬ 
tual  picked  SSR  times  (open  circles). 


value  T,.  The  GA  process  takes  1 030  s  for  one  shot  gather,  using  a 
1 .4-GHz  personal  computer.  The  typical  minimum  fitness  value  T,  is 
less  than  40  ms.  Thus,  the  average  SSR  residual  (7,/48)  between  the 
GA  inverted  time  and  the  picked  seismogram  time  is  less  than  1  ms, 
which  is  smaller  than  the  picking  error  and  represents  appropriate 
convergence  of  the  inversion  process. 

For  shot  2230  of  survey  line  DTIO,  the  GA  inversion  results  are 
compared  to  those  using  a  linear  cable  geometry  assumption.  A  more 
naturally  curving  cable  was  achieved  using  the  inverse  method  (sol¬ 
id  line.  Figure  4a).  Note  that  the  cable  depths  varied  from  850  m  at 
the  first  receiver  to  930  m  at  the  last  receiver  over  a  horizontal  dis¬ 
tance  of  approximately  420  m,  and  the  35-m  mismatch  between  the 
GA  inverted  geometry  and  the  linear  cable  assumption  produced  a 
34-ms  difference  in  corresponding  SSR  times  (Figure  4b).  The 
DTAGS  cable  was  closer  to  the  seafloor  than  to  the  sea  surface,  so  the 
time  difference  between  the  GA  method  and  linear  assumption  will 
be  even  larger  for  sediment  reflections.  Considering  that  a  1  -ms  time 
difference  will  produce  a  noticeable  mismatch  for  DTAGS  imaging, 
the  GA  inverted  cable  geometry  is  more  reliable  than  the  one  using  a 
linear  assumption. 


Datum  correction  and  NMO  calculation 

Once  the  cable  geometry  is  accurately  known  for  each  shot,  time 
corrections  are  computed  to  move  reflection  arrivals  on  the  original 
curved  array  to  the  corresponding  times  on  a  simulated  horizontal 
system.  We  assume  straight  raypaths  to  a  subseafloor  sediment  inter¬ 
face,  which  should  be  an  adequate  approximation.  Figure  5  illus¬ 
trates  that  the  original  raypath  AB  +  BC  from  the  source  to  the  last 
receiver  could  be  approximately  replaced  by  the  raypath  A'B' 
+  B'C',  with  both  source  and  receiver  on  a  horizontal  datum  at 
depth  //rf,  given  by  the  mean  value  of  source-receiver  depths,  and 
with  the  source-receiver  distance  corrected  to  the  original  offset  jc. 


Figure  5.  Illustration  of  time  correction  for  source/last  receiver  pair. 
The  original  raypath  AB  +  BC  from  source  to  the  last  receiver  re¬ 
flected  at  point  B  on  the  sediment  interface  is  replaced,  after  time 
correcting,  by  the  new  raypath  A'B'  +  B'C'  for  which  the  source 
and  receiver  are  at  the  same  datum  at  their  mean  depth  and  their  dis¬ 
tance  is  corrected  to  their  original  offset  on  the  cable.  V  is  the  aver¬ 
age  velocity  calculated  previously,  estimated  from  the  water  veloci¬ 
ty  and  a  sediment  velocity  structure  assumption;  Tonginai  is  the  origi¬ 
nal  twt  of  AB  -I-  BC;  7correc»ed  is  the  twt  of  A'B'  +  B'C'  after  time 
correction;  Ha  and  Hr  are  the  depth  of  the  datum  and  sediment  reflec¬ 
tor,  respectively,  measured  from  sea  surface;  dj  and  dr  are  the  depth 
of  the  source  and  the  receiver,  respectively. 
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The  raypath  changes  for  varying  incident/refiection  angles  at  dif¬ 
ferent  reflection  depths,  so  the  time  correction  is  calculated  by 
marching  a  small  time  window  iteratively  through  the  seismic  trace 
data.  The  previous  study  by  Walia  and  Hannay  (1999)  showed  that  a 
window  length  of  20  ms  is  adequate.  At  each  iteration,  the  depth  Hr 
of  the  reflector,  corresponding  to  the  two-way  time  Tongmai  for  ray- 
path  AB  +  BC,  is  calculated  assuming  an  approximate  velocity  pro¬ 
file  (e.g.,  from  nearby  multichannel  or  refraction  surveys).  This 
depth  Hr  is  then  used  to  predict  the  expected  time  Tcorrccwd  for  the  hor¬ 
izontal  source-receiver  system  at  the  datum  depth  H^i.  Time  correc¬ 
tions  (Tcorrectcd)  for  the  top  and  bottom  of  the  marching  window  are 
calculated  separately,  and  the  data  in  the  original  window  are 
stretched/contracted  and  time  shifted  to  fit  the  new  window  accord¬ 
ingly.  The  time  correction  calculations  are  differential  in  nature,  so 
the  errors  are  only  weakly  dependent  on  the  velocity  profile  (Walia 
and  Hannay,  1999),  particularly  because  the  DTAGS  signal  penetra¬ 
tion  is  limited  to  the  upper  400  m  of  sediments.  Thus,  it  is  adequate 
to  use  a  very  approximate  velocity  profile,  such  as  a  constant  average 


velocity  over  the  raypaths  estimated  by  combining  the  water  and 
shallow  sediments  (e.g.,  1600  m/s). 

Unlike  the  constant  datum  depth  used  for  the  whole  survey  line  by 
Walia  and  Hannay  (1999),  an  individual  datum  was  chosen  at  the 
mean  depth  of  each  source  and  receiver  pair.  That  is,  a  floating  datum 
system  was  used  in  this  study  not  only  for  different  shot  instants  but 
also  for  different  receivers  within  the  same  shot  gather. 

Once  the  time  corrections  for  source-receiver  pairs  have  been  ap¬ 
plied,  an  NMO  process  is  used  to  shift  traces  to  equivalent  vertical 
incidence  using  constant  seawater  velocity.  For  the  DTAGS  imag¬ 
ing,  we  mainly  used  the  nearer-offset  dense  acoustic  array  for  which 
the  moveout  is  quite  small  and  using  water  velocity  provides  a  rea¬ 
sonable  image,  especially  for  seafloor  alignment.  For  comparison 
with  reflection  images  collected  with  surface  source  and  receivers, 
the  DTAGS  data  are  finally  corrected  to  put  shot  and  receiver  at  the 
sea  surface  by  shifting  the  trace  downward  by  a  two-way  time  at  sea¬ 
water  velocity  for  the  vertical  raypath  from  the  sea  surface  to  the  da¬ 
tum  depth.  A  few  shot  gathers  are  shown  in  Figure  6  before  and  after 


Figure  6.  Shot  gathers  from  the  deep-tow  array  be¬ 
fore  (top)  and  after  (bottom)  time  corrections  for 
variations  in  both  receiver  array  shape  and  shot 
depth.  The  bottom  images  have  been  corrected  to 
the  vertical  traveltime  from  the  sea  surface  to  the 
seafloor,  and  the  traces  are  not  equally  spaced  in 
offset:  the  first  24  channels  are  the  acoustic  array 
(channel  124,  spacing  3  m),  and  the  remaining  24 
channels  with  much  larger  moveout  corresponds  to 
the  geophysical  array  with  longer  offset  (channel 
2548,  spacing  15  m)  (see  Table  1 ). 


Shot  256  Shot  257  Shot  258 


Downloaded  21  Jul  2009  to  128.160.51.126.  Redistribution  subject  to  SEG  license  or  copyright;  see  Terms  of  Use  at  http://segdl.org/ 


Imaging  a  hydrate-related  cold  vent 


B29 


time  corrections  for  variations  in  receiver  array  shape  and  source 
depth;  the  true  reflection  hyperbolas  have  been  restored,  and  the 
sources  are  at  a  consistent  depth. 

DTAGS  data  imaging 

To  form  the  high-resolution  seismic  image,  the  time-corrected 
DTAGS  data  were  passed  to  a  receiver  statics  analysis  and  a  shot  im¬ 
age  stitching  process.  To  seek  consistency  along  a  seismic  reflector, 
these  calculations  were  based  on  crosscorrelations  of  trace  enve¬ 
lopes. 

To  remove  surface-consistent  residual  statics,  receiver  statics 
were  first  applied  to  minimize  the  difference  between  adjacent  re¬ 
ceivers.  Fora  given  shot,  the  trace  envelopes  of  each  adjacent  receiv¬ 
er  pair  were  crosscorrelated  and  the  lag  between  the  traces,  corre¬ 
sponding  to  the  maximum  amplitude  of  the  crosscorrelation,  was  de¬ 
termined.  The  mean  value  of  these  lags  over  all  shots,  after  conver¬ 
sion  to  sample  number  or  static  time,  was  the  relative  statics  correc¬ 
tion  for  these  two  receivers.  The  calculated  receiver  statics  were  very 
small  (<  1  ms)  because  the  maximum  15-m  offset  between  receiv¬ 
ers  is  small  compared  to  the  source  Fresnel  zone. 

The  major  differences  between  two  adjacent  shot  images  are 
caused  mainly  by  variations  in  the  source  spacing  and  source  depth. 
Although  the  navigation  of  the  ship  provided  estimates  of  ship  posi¬ 
tion  and  speed,  the  source  did  not  follow  the  motion  of  the  ship  pre¬ 
cisely.  With  the  source  located  at  the  end  of  a  long  tow  cable,  the 
source  would  rise  or  fall  in  the  water  column  as  the  speed  of  the  ship 
varied,  and  the  horizontal  distance  of  the  source  behind  the  ship  also 
varied.  To  constrain  this  variation  in  source  position,  the  large  over¬ 
lapping  area  between  two  adjacent  shot  images  was  used  to  calculate 
vertical  shot  statics  and  the  optimal  horizontal  position  for  image 
stitching. 

Although  the  nominal  ship  speed  was  approximately  2  knots, 
such  low  speeds  are  difficult  to  maintain  at  a  con¬ 
stant.  The  ship  speed  varied  from  1 .5  to  3  knots 
over  most  of  the  survey,  with  corresponding  shot 
spacing  variations  of  2145  m  for  a  nominal  shot 
interval  of  30  s.  Because  the  maximum  trace  off¬ 
set  is  438  m  (Table  1),  the  reflection  midpoints 
extend  out  to  219  m  from  the  source.  There  is 
thus  a  significant  overlap  in  common-reflection- 
point  (CRP)  position  from  one  shot  to  the  next, 
varying  from  174  to  198  m. 

The  overlapping  areas  for  adjacent  shot  images 
were  compared  CRP  by  CRP  for  shot  spacing 
varying  from  21  to  45  m  with  a  stepping  size  of 
3  m  (corresponding  to  the  spacing  of  the  acoustic 
array).  This  produced  a  series  of  crosscorrelation 
seismograms  for  different  overlapping  areas, 
with  each  series  characterized  by  a  mean  value 
and  a  standard  deviation  over  all  the  picked  trace 
lags:  the  sample  numbers  corresponding  to  the 
maximum  amplitude  of  the  crosscorrelation.  The 
optimal  shot  spacing,  providing  the  most  simi¬ 
lar  overlapping  image  between  adjacent  shots,  is 
the  spacing  with  the  smallest  mean  lag  value. 

In  calculating  the  mean  lag,  we  rejected  shot 
spacings  for  which  the  standard  deviation  was 
more  than  twice  the  smallest  value  of  all  standard 
deviations. 


A  100%  coverage  image  is  formed  using  selected  traces  from  the 
acoustic  array  only.  The  farther  offset  geophysical  array  traces  and 
overlapping  acoustic  array  traces  are  used  for  velocity  analysis  only 
and  are  not  stacked  in  the  image.  For  these  larger-offset  traces,  errors 
in  alignment  and  in  NMO  correction  are  larger  than  for  the  near-off¬ 
set  acoustic  array,  so  the  single-fold  image  has  the  greatest  continu¬ 
ity  in  image  stitching  and  the  highest  resolution  in  image  quality. 

The  above  process,  using  crosscorrelation  lags  among  all  CRP 
traces  between  adjacent  shots,  defines  the  optimal  overlap  in  traces 
or  the  optimal  horizontal  shot  position.  To  calculate  the  relative  ver¬ 
tical  shot  statics  correction  for  the  1 00%  coverage  image,  only  lags 
for  the  first  three  common  traces  were  used.  This  process  automati¬ 
cally  stitched  80%-90%  of  the  DTAGS  shot  images.  Some  addition¬ 
al  manual  adjustment  for  individual  shots  was  required,  particularly 
for  low-grade  shots  where  many  empty  or  noisy  traces  existed. 

There  are  two  obvious  advantages  of  this  image  stitching  method. 
First,  calculating  relative  shot  spacing  does  not  rely  on  any  position¬ 
ing  information  (e.g.,  GPS  navigation,  acoustic  transponder  net). 
Second,  a  more  reliable  CRP  grouping  could  be  achieved  for  the 
NMO  velocity  analysis. 

The  stitched  DTAGS  images  using  the  above  method  have  no  in¬ 
formation  on  geographic  coordinates.  To  better  correlate  subseafloor 
features  to  seafloor  topography,  a  stitched  DTAGS  image  was 
aligned  to  the  correspx)nding  rover  track  by  assigning  the  beginning 
and  end  of  the  image  to  the  beginning  and  end  of  the  rover  track,  and 
interpolating  trace  coordinates  between  the  nearest  rover  coordi¬ 
nates.  We  estimate  errors  in  absolute  positions  are  about  2030  m, 
based  on  comparison  with  a  single-channel  seismic  (SCS)  grid  col¬ 
lected  in  2000.  Figure  7  is  an  example  of  the  final  processed  acoustic 
image  of  the  2002  DTAGS  survey  line  DT09. 
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Figure  7.  Acoustic  image  (100%  coverage)  of  the  survey  line  DT09,  showing  high  resolu¬ 
tion  of  subseafloor  structures  laterally  and  vertically.  The  location  of  lODP  drill  site 
U 1 328C  is  labeled.  S3 1 7  and  S377  label  two  locations  where  conventional  velocity  anal¬ 
ysis  was  carried  out. 
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Moveout  velocity  analysis 

In  the  first  attempt  at  velocity  analysis  of  DTAGS  data,  Rowe  and 
Gcttrust  (1993)  used  shot  gathers.  Velocity  uncertainties  were  quite 
high  because  reflection  points  extended  over  a  distance  of 
>200  m. Sorting  by  common  reflection  point  provides  significant 
improvement  in  precision,  but  difficulties  arise  because  of  the  irreg¬ 
ular  trace  spacing  resulting  from  variable  ship  speed  and  array  ge¬ 
ometry.  In  the  DTAGS  velocity  study  by  Chapman  et  al.  (2002),  the 
shot  spacing  was  estimated  approximately  only  for  line  segments 
where  ship  speed  was  relatively  consistent. 

In  the  current  study,  more  accurate  CRP  sorting  was  achieved  by 
better  determination  of  the  optimal  shot  spacing  from  stitching  adja¬ 
cent  shot  images.  The  CRP  bin  size  was  set  to  7.5  m,  corresponding 
to  the  15-m  group  spacing  in  the  geophysical  array.  However,  to  in¬ 
crease  signal-noise  ratio  for  velocity  analysis,  four  adjacent  CRPs 
were  grouped  to  form  a  CRP  supergather.  The  effective  bin  size  was 
then  '-SO  m,  corresponding  to  the  nominal  shot  spacing.  So  there  are 
as  many  offsets  in  a  CRPgather  as  in  a  shot  gather. 

Similar  to  the  image  process,  time  corrections  for  cable  geometry 
and  source-depth  variations  are  required  to  obtain  a  horizontal  sys¬ 


tem.  For  each  CRPsupergather,  however,  only  a  single  reference  da¬ 
tum  is  chosen  at  the  average  value  of  mean  depths  of  each  source-re¬ 
ceiver  pair;  this  is  adequate  for  velocity  analysis  because  the  source 
depth  varies  by  only  a  small  amount  over  several  successive  shots. 
No  depth  shifting  to  the  sea  surface  is  needed.  Once  the  cable  geome¬ 
try  is  corrected,  the  CRP  supergather  can  be  treated  approximately  as 
a  conventional  seismic  data  set,  and  standard  seismic  analysis  rou¬ 
tines  can  be  used. 

There  are  inadequate  reflector  horizons  within  the  main  vent 
structure  for  velocity  determination,  and  velocity-depth  estimates 
were  obtained  from  CRP  supergathers  S3 17  and  S337  of  DTAGS2 
survey  line  DT09  (Figure  7  and  Figure  13)  because  they  contain 
strong  reflections  that  are  nearly  horizontal.  They  are  centered  ap¬ 
proximately  1 80-270  m  to  the  southwest  and  northeast  sides  of  the 
center  of  the  prominent  seafloor  pxxkmark  (distance  ~2460  m)  as¬ 
sociated  with  Bullseye  vent.  The  reflection  points  for  the  CRP  super¬ 
gathers  span  about  25  m,  from  2175  to  2200  m  for  S3 17  and  from 
2620  to  2645  m  for  S337. 

Figure  8a  shows  prestack  images  of  the  constant-velocity  gathers 
for  CRP  supergather  S3 1 7.  For  a  given  reflection,  we  picked  the  rms 
velocity  that  best  flattens  the  reflection.  After  applying  NMO  to  the 
data,  the  final  prestack  image  in  the  last  panel 
shows  all  flattened  reflections.  From  the  mono- 
tonically  increasing  rms  velocity  profile,  the  cor¬ 
responding  interval  velocities  can  be  calculated 
using  the  Dix  equation  (Figure  9a).  The  deepest 
reflection  occurs  at  195  mbsf,  slightly  shallower 
than  the  expected  depth  of  2 1 0  mbsf  for  the  BSR. 
For  CRPsupergather  S377,  similar  NMO  veloci¬ 
ty  panels  were  generated  and  interval  velocities 
calculated  (Figures  8b  and  9b).  The  velocity  anal¬ 
ysis  was  limited  to  the  upper  125  m  because  of 
weak  reflectors  at  greater  depth. 

Velocity  uncertainty  estimation 

The  uncertainties  of  rms  velocities  from  the 
constant-velocity  gathers  were  about  ±2.5  m/s 
given  the  velocity  step  size  of  5  m/s.  The  associ¬ 
ated  uncertainties  for  interval  velocities  (thick 
solid  lines  in  Figure  9)  are  significantly  larger 
(±50  m/s  to  >  100  m/s),  particularly  for  thin 
sediment  layers.  The  errors  are  in  part  caused  by 
limitations  in  the  datum  correction  and  in  the  hy¬ 
perbolic  assumption  and  the  associated  assump¬ 
tion  that  interval  velocities  can  be  calculated  from 
rms  velocities  using  the  Dix  equation.  During  the 
datum  correction  (Figure  5),  the  true  raypath 
ABC  was  replaced  by  an  approximate  one 
A'B'C',  and  this  approximation  produced  errors 
in  source-receiver  offsets  and  vertical  time  shifts. 
The  largest  errors  occurred  for  the  near-seafloor 
sediments,  and  the  error  decreased  quickly  as  the 
depth  of  reflector  increased.  Other  contributions 
to  error  are  limitations  in  parameter  estimation, 
notably  source/receiver  geometry  and  correction 
velocities  for  the  water  layer  and  sediment  layers. 

For  reflection  alignment  in  velocity  analysis  or 
for  reflection  stacking,  the  depth  differences  be¬ 
tween  the  source  and  receivers  should  be  less  than 
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Figure  8.  Constantvelocity  gathers  (CVG)  at  CRPsupergather  (a)  S31 7  and  (b)  S377  on 
DTAGS2  survey  line  DT09.  The  rms  velocity  (in  m/s)  is  shown  at  the  top  of  each  panel, 
and  the  bold  arrow  indicates  the  specific  event  that  is  best  flattened  with  that  velocity.  The 
location  where  the  acoustic  array  changes  to  the  geophysical  array  is  also  indicated  by  a 
vertical  arrow  at  the  top  of  the  first  panel.  The  ALL  panels  show  the  resulting  images  after 
NMO  with  the  final  rms  velocity  profiles  (Figure  9). 
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DTAGS  OBSERVATIONS  AT 
BULLSEYE  VENT 

Imaging  of  hydrate  cap 

In  July  2000,  piston  cores  recovered  massive 
hydrate  at  depths  of  38  m  below  the  seafloor 
(small  dark  gray  patch  within  the  light  gray  Bulls- 
eye  vent  area;  see  Figure  2).  On  a  grid  of  SCS 
lines,  a  shallow  reflector  with  a  cap  shape  was  in¬ 
terpreted  as  the  top  of  massive  hydrate  on 
crossline  XL75,  oriented  southeast-northwest 
(Figure  1 0;  see  Figure  2  for  location).  The  hydrate 
cap  deepens  away  from  the  center  of  the  vent,  and 
its  depth  matches  well  with  the  depths  of  the  hy¬ 
drate  recovery  at  several  piston  core  sites  (Riedel, 
2001 ;  Riedel  et  al.,  2002, 2006b). 

The  hydrate  cap  reflector  is  also  seen  on  the 
2002  DTAGS  lines.  Lines  DTI  6  and  DTI 7  cross 
the  Bullseye  vent  and  are  oriented  approximately 
southeast-northwest,  close  to  the  direction  of  the 
SCS  crossline  XL75  but  offset  by  50-100  m  (lo¬ 
cations  on  Figure  2).  Also,  DTAGS  lines  DT08, 
DT09,  and  DTIO  cross  the  center  of  the  vent 
southwest-northeast  (Figure  2).  On  these  lines, 
the  top  of  the  hydrate  cap  (black  curves.  Figures 
11-13)  was  inferred  from  individual  images  and 
from  intersections  with  other  DTAGS  lines.  For 
comparison,  the  top  of  the  hydrate  cap  was  also 
inferred  from  the  SCS  grid  (green  curves).  The 
DTAGS  data  and  the  SCS  data  do  not  give  exactly 
the  same  depths  for  the  hydrate  cap;  this  is  be¬ 
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Figure  9.  The  rms  velocity  (Yms)  and  corresponding  interval  velocity  (Vj  in  two-way 
traveltime  and  depth  below  seafloor,  determined  from  constant-velocity  gathers  forCRP 
supergather  (a)  S3 1 7  and  (b)  S77  in  Figure  8.  The  estimated  error  bounds  for  Vint  are  indi¬ 
cated  by  light  gray  shading.  The  no-hydrate  reference  velocity  of  Chen  (2006)  is  marked 
by  the  dashed  line. 


half  the  wavelength  of  the  signal.  Thus,  for  the  high-resolution 
DTAGS  with  a  midfrequency  of  500  Hz,  the  vertical  positions  of  the 
source  and  receivers  must  be  known  to  within  1 .5  m.  The  errors  in 
the  array  geometry  correction  are  difficult  to  quantify,  but  we  can  es¬ 
timate  the  sensitivity  of  the  velocity  analysis  to  variations  in  source- 
receiver  offset  from  the  Norwegian-Greenland  Sea. 

With  the  inverted  cable  geometry  of  shot  196  on  line  DT09,  we 
created  a  simple  1 D  model  with  water  depth  1 270  m,  water  velocity 
1480  m/s,  and  uniform  sediment  velocity  1550  m/s.  Using  SnelVs 
law,  we  traced  rays  from  the  source  to  the  receivers  from  a  sediment 
interface  at  100  mbsf.  With  these  synthetic  traveltimes  as  input,  we 
applied  our  datum  correction  followed  by  our  NMO  velocity  analy¬ 
sis  to  calculate  the  sediment  interval  velocity.  We  then  varied  an  indi¬ 
vidual  model  parameter  to  determine  the  effect  on  the  inverted  sedi¬ 
ment  velocity.  We  found  a  5-m  error  in  the  maximum  source-receiv¬ 
er  offset  caused  a  large  interval  velocity  change  of  20  m/s.  We  also 
discovered  a  10-m/s  uncertainty  in  water  velocity  only  causes  an  in¬ 
terval  velocity  change  of  <  3  m/s.  Another  discovery  was  that  an  un¬ 
certainty  of  50  m/s  in  the  sediment  velocity  used  for  datum  correc¬ 
tion  only  causes  a  change  of  1 0  m/s  in  the  invert¬ 
ed  interval  velocity.  Finally,  we  learned  a  very 
large  bias  of  100  m  from  zero-error  datum  only 
causes  an  interval  velocity  change  of  <20  m/s. 

Thus,  the  most  sensitive  parameters  for  determin¬ 
ing  sediment  velocity  are  those  associated  with 
cable  geometry  and  datum  correction. 


cause  their  dominant  frequencies  are  significantly  different  (400  Hz 
for  DTAGS,  1 20  Hz  for  SCS),  so  amplitudes  vary  from  different  at¬ 
tenuation  and  different  respx^nses  to  velocity  gradients  and  thin  lay¬ 
ers. 

lODP  site  U1328  was  drilled  very  close  to  line  DT09  (Figure  2) 
within  the  Bullseye  vent.  DT09  shows  highly  resolved  details  of  the 
sediment  structure  (Figure  1 3).  The  main  vent  is  indicated  by  the 
seismic  wipeout  zone  between  2250-2350  m  distance.  This  blank 
zone  disrupts  the  set  of  strong  coherent  reflectors  that  occur  at 
~1 525  ms  below  the  seafloor;  however,  it  does  not  disrupt  reflectors 
at  the  seafloor.  There  is  a  high  reflectivity  zone  dipping  southwest 
from  a  seabed  depression,  where  the  hydrate  cap  nearly  intersects 
the  seafloor,  down  to  a  depth  of  2025  ms  two-way  time  (twt) 
('^1025  mbsf)  where  it  meets  the  blank  zone.  This  high  reflectivity 
feature  may  be  a  fracture  zone  that  provides  an  oblique  pathway  for 
fluid  flow  in  the  main  vent  of  the  blank  zone  to  bring  gas  and/or  hy¬ 
drate-coated  bubbles  to  very  near  the  seafloor.  The  top  of  the  fracture 
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Line  DTI  6  continues  as  a  long  transit  line  extending  to  the  north¬ 
west.  The  1999  COAMS  (Canadian  Ocean  Acoustic  Measurement 
System)  grid  inline  IN26  is  coincident  with  DTI  6  (Figure  1 ).  A  com¬ 
bination  of  the  surface-towed  seismic  data  and  the  deep-towed 
DTAGS  data  is  needed  to  provide  sufficient  information  for  propter 
imaging  and  interpretation.  The  high-quality  seismic  image  of  DTI  6 
shows  very  detailed  shallow  sediment  layers  (Figure  14b),  com¬ 
pared  to  the  conventional  multichannel  seismic  image  of  IN26  (Fig¬ 
ure  14a).  The  1N26  image  resolves  four  seismic  wipeout  zones  as  la¬ 
beled  in  the  slope  sediments,  likely  representing  fluid  vent  sites  or 
proto  vents  because  not  all  of  them  extend  upward  to  the  seafloor. 
However,  on  the  DTI 6  image,  small  diffractions  are  present  at  the 
seafloor  for  blank  zones  BZ2,  BZ3,  and  BZ4,  indicating  the  surface 
expressions  of  these  ventlike  features.  The  surface  expression  of  the 
largest  blank  zone  BZl  (Bullseye  vent)  on  DTI  6  is  similar  to  a  pock¬ 
mark  with  a  depression  in  the  middle  (Figure  14  and  Figure  1 3b). 

The  second  set  of  strong  arrivals  at  about  20  ms  two- way -time  be¬ 
low  the  seafloor  seems  an  apparent  upward  barrier  for  blank  zone 
BZ2  (Figure  13b);  the  small  pull-up  or  deformation  that  indicates  the 
possible  vent  is  restricted  to  depths  below  the  strong  reflectors. 
However,  blank  zones  BZ3  and  BZ4  partially  cross  this  layer.  The 
flanks  of  the  blank  zones  show  small-scale  diffractions  and,  for 
blank  zone  BZ3,  an  apparent  increase  in  seismic  reflectivity  for  sev¬ 
eral  reflections  just  outside  the  vent  (Figure  14b). 

The  entire  region  from  BZ2  to  just  northwest  of 
BZ4  (Figure  14b)  is  a  structural  high.  At  the  high- 
31TO  ^  3200  frequency  range  of  DTAGS,  the  amplitude  pat¬ 

tern  in  this  region  also  contrasts  with  the  sur¬ 
rounding  area:  reflections  from  about  1.8-1.87  s 
twt  have  enhanced  amplitude,  but  the  deepest 
prominent  reflection  between  BZl  and  BZ2 
(-^1.95  s  twt)  does  not  continue  beneath  the 
structural  high.  The  entire  structural  high  region 
may  contain  greater  amounts  of  hydrate  than  the 
surrounding  area.  Such  an  interpretation  is  con¬ 
sistent  with  the  observations  from  controlled- 
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Figure  10.  Southeast-northwest  seismic  section,  showing  hydrate 
cap  details  from  Teledyne  2000  crossline  XL75  (see  Figure  2  for  lo¬ 
cation).  The  locations  of  intersection  with  DTAGS  line  DT09  and 
nearby  lODP  site  U1 328  are  also  labeled.  Note  reduced  seismic  am¬ 
plitudes  below  the  hydrate  cap  (center  of  the  vent). 


zone  (black  curve.  Figure  1 3)  corresponds  to  the  top  of  the  hydrate 
cap. 

Blank  zones  on  transit  line  DT16 


Figure  1 1 .  Southeast-northwest  seismic  section,  showing  hydrate  cap  details  from  DTI  6. 
The  green  curve  is  the  hydrate  cap  location  indicated  from  the  Teledyne  2000  grid;  the 
black  curve  is  from  the  DTAGS  line.  The  red  vertical  line  indicates  the  intersection  loca¬ 
tion  with  line  DT09. 


source  EM  measurements  (Schwalenberg  et  al, 
2005),  in  which  electrical  resistivity  beneath  the 
entire  BZ2  to  BZ4  region  is  higher  than  that  in  the 
no-hydrate  reference  area. 

A  prominent  small  ridge  lies  at  the  foot  of  the 
slope,  approximately  1 500  m  northwest  of  blank¬ 
ing  zone  BZ4.  Significant  seismic  blanking  oc¬ 
curs  beneath  the  ridge.  On  neighboring  COAMS 
high-frequency  surface  multichannel  lines,  the 
ridge  with  its  associated  blanking  is  buried  at 
greater  depths  beneath  less-deformed  slope  sedi¬ 
ments  (Riedel,  2001).  Also,  the  DTAGS  image 
DTI 6  shows  a  clear  second  coherent  interface 
and  a  diapiric  structure  with  uplifted  sediments 
on  both  sides  below  the  ridge.  Thus,  the  blanking 
beneath  the  ridge  may  be  the  result  of  widespread 
penetrative  deformation  in  the  accreted  sedi¬ 
ments  (Davis  and  Hyndman,  1989)  and  not  from 
gas-charged  or  hydrate-filled  sediments. 
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Figure  12.  Southeast-northwest  seismic  section, 
showing  hydrate  cap  details  from  DTI  7.  The  green 
curve  is  the  hydrate  cap  location  indicated  from  the 
Teledyne  2000  grid;  the  black  curve  is  from  the 
DTAGS  line.  The  red  vertical  line  indicates  the  in¬ 
tersection  location  with  line  DT09. 
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Figure  13.  Southwest-northeast  seismic  section 
from  DTAGS2  line  DT09.  The  green  curve  is  the 
hydrate  cap  location  indicated  from  the  Teledyne 
2000  grid;  the  black  curve  is  from  the  DTAGS  line. 
The  dashed  black  rectangle  box  indicates  a  layer 
with  high  hydrate  concentrations  from  about 
540  mbsf  revealed  by  the  resistivity  log  of  lODP 
site  1 328  hole  A  on  this  line  (Figure  1 5).  The  top  of 
this  layer  corresponds  to  the  hydrate  cap  location. 
Maximum  hydrate  concentrations  ('^60%)  corre¬ 
spond  to  a  high-reflectivity  fracture  zone  that 
projects  to  the  seafloor.  The  two  locations  where 
velocity  analyses  were  applied  are  also  labeled 
(S317,S377). 
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Velocity  results  and  comparison  with  lODP  drilling 

One  of  the  objectives  of  the  2002  DTAGS  survey  was  to  investi¬ 
gate  the  properties  of  gas  hydrate  and  associated  BSRs.  The  P-wave 
seismic  velocity  of  pure  methane  hydrate  is  3650  m/s  (Waite  et  al., 
2000),  significantly  higher  than  the  typical  1500-1600-m/s  seismic 
velocity  of  sediments  in  the  upper  200-300  mbsf  if  no  gas  hydrate  is 
present.  Thus,  the  measured  seismic  velocity  of  sediments  contain¬ 
ing  hydrate  provides  a  quantitative  measure  of  the  concentration  of 
hydrate  in  the  fluid  pore  space.  Because  of  its  deep-towed  geometry 
and  high-frequency  content,  DTAGS  is  well  suited  in  principle  to 
provide  high -resolution,  well-constrained  velocity  information 
about  sediments  containing  hydrate. 

At  the  two  test  supergather  locations  S317  and  S377  (Figure  15), 
NMO-derived  interval  velocities  can  be  compared  to  the  sonic  log 
and  vertical  seismic  profile  (VSP)  from  lODP  site  U1328  located 
within  the  Bullseye  vent.  The  supergathers  are  located  just  outside 
the  seismic  wipeout  region  that  characterizes  the  vent,  less  than 
250  m  from  U1328  (Figure  7);  unfortunately,  because  of  the  blank¬ 
ing,  NMO  velocities  cannot  be  determined  from  within  the  vent  it¬ 
self  because  there  are  few  clear  reflectors.  Also,  downhole  velocity 
information  is  available  at  site  U1328  only  for  depths  greater  than 
about  60  m.  NMO  velocities  have  been  measured  at  locations  just 
outside  the  vent,  from  shallow  depths  (<60  mbsf)  where  no  sonic 
velocities  are  available,  down  to  depths  that  approach  the  BSR. 


At  supergather  S3 1 7  just  southwest  of  the  Bullseye  vent,  interval 
velocities  increase  from  1 .53  km/s  in  the  first  30-m-thick  layer  be¬ 
low  the  seafloor  to  1.67  km/s  for  the  deepest  layer  that  extends  to 
195  mbsf  (Figure  15).  These  velocities  are  very  close  to  the  conven¬ 
tional  multichannel  estimate  by  Chen  (2006)  of  the  reference  seis¬ 
mic  velocity  depth  in  this  region  for  sediments  containing  no  hy¬ 
drate.  The  velocities  are  generally  consistent  with  the  sonic  veloci¬ 
ties  from  site  U1328  and  with  the  VSP  values  measured  at  depths 
>  100  mbsf,  which  indicate  an  average  velocity  of  1645  m/s  (Fig¬ 
ure  1 5).  In  the  two  layers  extending  from  055  mbsf,  there  are  normal 
velocities  with  no  indication  of  high-velocity  hydrate.  At  nearby  site 
U1328,  very  high  logging- while-drilling  (LWD)  resistivities  were 
measured  at  depths  of  540  mbsf  (Figure  15),  interpreted  as  a  region 
of  massive  hydrate  with  saturations  of  60%-80%  of  the  pore  space 
(Riedelet  al.,  2006a;  Chen,  2006).  This  velocity  analysis  indicates 
that  the  shallow  massive  hydrate  does  not  extend  outside  the  vent  to 
the  southwest. 

At  supergather  S377  on  the  northeast  side  of  the  seafloor  pock¬ 
mark,  the  derived  interval  velocities  for  most  layers  are  consistent, 
within  the  estimated  uncertainties,  with  the  no-hydrate  reference  ve¬ 
locity  profile  of  Chen  (2006)  and  with  the  measured  sonic  velocities 
at  site  U1 328  below  a  depth  of  about  60  mbsf  (Figure  1 5).  However, 
the  interval  velocity  from  50-70  mbsf  is  1.69  km/s.  Although  the 
estimated  error  bounds  extend  to  just  include  the  reference  profile 
value,  the  preferred  value  implies  a  possible  large  velocity  increase 
at  this  depth,  with  velocities  1 00  m/s  greater  than  the  reference.  Us- 
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Figure  1 4.  (a)  Conventional  multichannel  seismic 
image  of  the  1 999  COAMS  grid  inline  IN26,  pro¬ 
viding  a  good  comparison  to  the  coincident 
DTAGS  transit  line  DTI  6.  (b)  Acoustic  image 
(100%  coverage)  of  the  2002  DTAGS2  transit  line 
DTI 6,  showing  fine  layering  of  slope  sediments, 
four  seismic  blanking  zones  BZl-4,  and  a  small 
ridge  (SR)  of  accreted  sediments,  possibly  diapiric. 
The  dashed  line  indicates  the  boundary  between  the 
less-deformed  region  (LDR)  of  slope  sediments 
and  the  highly  deformed  region  (HDR)  of  older  ac¬ 
creted  sediments. 
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Figure  15.  Comparisons  of  moveoul  velocity  analysis  results  with 
U1328  drillhole  logs.  Velocity  panel:  sonic  log  (blue  line),  interval 
velocities  of  super-CRP  S3 17  (thick  green  line)  and  S377  (thick  red 
dashed  line)  (error  bounds  given  in  Figure  9),  Bayesian  inversion  re¬ 
sults  of  VSP  starting  from  1 06  mbsf  (thick  black  line;  associated  er¬ 
rors  usually  a  minimum  of  ±2{X)  m/s),  and  the  no-hydrate  reference 
velocity  of  Chen  (2(X)6)  (thick  black  dashed  line).  Resistivity  panel: 
Button  deep  resistivity  average  (BDAV,  red  line)  and  phase  shift  re¬ 
sistivity  blended,  1 6  inxhwa  (P16B,  blue  line).  TTie  sonic  log  was  re¬ 
corded  in  hole  U1328C,  whereas  resistivity  logs  came  from  nearby 
drillhole  U 1328 A,  38  m  northeast  of  U 1 328C  (Riedel  etal.,  2006a). 

ing  the  velocity-concentration  models  developed  by  Chen  (2{X)6)  for 
nearby  site  U 1 327,  we  determine  that  the  average  hydrate  concentra¬ 
tion  in  this  layer  is  1 5%-20%. 

CONCLUSIONS 

The  deep- towed  DTAGS  system  has  proven  capable  of  imaging 
shallow  structure  beneath  the  seafloor  at  much  higher  resolution 
than  surface-towed  multichannel  systems,  but  careful  data  analysis 
is  required.  A  data  analysis  process  was  developed  to  determine  the 
geometry  of  the  receiver  cable  for  each  shot  and  to  determine  the  po¬ 
sition  of  the  shots  along  a  profile.  Cable  geometry  is  determined  by 
inverting  known  node  depths  and  picked  traveltimes  for  the  direct 
arrivals  and  reflections  from  the  sea  surface.  For  source  positions, 
envelopes  of  overlapping  traces  between  two  adjacent  gathers  are 
cross  correlated  to  find  the  optimal  spacing  between  the  shots.  Al¬ 
though  all  traces  along  the  438 -m  cable  are  used  for  the  source-re¬ 
ceiver  geometry  and  for  subsequent  NMO-velocity  analysis,  the 
highest  resolution  structural  image  is  found  using  only  the  short-off- 
set  (78-m)  high-density  acoustic  array.  At  the  high-frequency  short 


wavelengths  of  the  DTAGS  system,  positioning  is  still  inadequate 
for  good  trace-to-trace  coherence  using  the  full  array. 

The  DTAGS  images  resolve  multiple  seismic  wipeout  zones  in 
the  midslope  region  of  the  Northern  Cascadia  margin.  The  blank 
zones  likely  represent  fluid  vent  sites.  A  few  may  be  protovents  be¬ 
cause  not  all  of  them  extend  upward  to  the  seafloor.  At  the  largest 
vent  in  the  region,  the  Bullseye  vent,  the  DTAGS  data  image  is  a  10 
-m- thick  high-reflectivity  zone,  extending  from  the  seafloor  to  a 
depth  of  about  30  m  and  dipping  seaward  (southwest)  at  a  small  an¬ 
gle  of  7®.  This  zone  likely  corresponds  to  the  upper  portion  of  the 
massive  methane  hydrate  layer  (~35  m  thick)  that  has  been  piston 
cored  and  encountered  at  lODP  site  U1328.  This  high-reflectivity 
zone  is  interpreted  as  a  complex  near-surface  hydraulic  fracture  sys¬ 
tem  filled  by  hydrate,  with  gas-rich  fluids  fed  from  steeply  dipping 
faults  below  it. 

Interval  velocities  from  NMO  analysis  of  datum-corrected 
DTAGS  data  indicate  that  the  shallow  layer  of  massive  hydrate  does 
not  extend  outside  the  vent  region  because  velocities  in  the  upper 
50  m  just  outside  the  vent  generally  are  consistent  with  the  reference 
velocity  profile  for  sediments  containing  no  hydrate.  However,  at  a 
depth  of  50-70  mbsf  just  to  the  northeast  of  the  vent,  interval  veloci¬ 
ties  are  approximately  100  m/s  higher  than  the  no-hydrate  refer¬ 
ence.  This  high-velocity  layer  may  be  a  fluid-channeling  fault  zone. 
The  high  velocities  suggest  that  hydrate  concentrations  in  this  layer 
are  1 5%-20%  of  the  pore  space. 
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